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Spherical dust collapse generally forms a shell focusing naked singularity at the symmetric center. 
This naked singularity is massless. Further the Newtonian gravitational potential and speed of the 
dust fluid elements are everywhere much smaller than unity until the central shell focusing naked 
singularity formation if an appropriate initial condition is set up. Although such a situation is 
highly relativistic, the analysis by the Newtonian approximation scheme is available even in the 
' vicinity of the space-time singularity. This remarkable feature makes the analysis of such singularity 

formation very easy. We investigate non-spherical even-parity matter perturbations in this scheme 
• by complementary using numerical and semi-analytical approaches, and estimate linear gravitational 

waves generated in the neighborhood of the naked singularity by the quadrupole formula. The result 
shows good agreement with the relativistic perturbation analysis recently performed by Iguchi et 
al. The energy flux of the gravitational waves is finite but the space-time curvature carried by them 
, diverges. 
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C-i '. PACS number(s): 04.25.Nx,04.30.Db,04.20.Dw 

I. INTRODUCTION 

X 

S-H 

General relativity predicts that gravitational collapse of massive objects leads to space-time singularities in rather 
general circumstances in our universe [jl] — ^| and such singularities might be accompanied with the blow up of physical 
quantities (energy density, pressure and curvature of space-time, etc). The known physical laws including general 
relativity itself will break down in the neighborhood of the space-time singularity and hence the quantum theory of 
gravity is believed to be necessary to describe the physical phenomena in such a region. 

One of the important issues is whether the space-time singularities formed in our universe are visible or not for the 
observer (us) far from the region where the gravitational collapse occurs. The cosmic censorship conjecture proposed 
by Penrose gave a strong motivation to investigate this problem j|. Roughly speaking, this conjecture states that the 
singularity is not visible for any observer if it is resulted from physically reasonable initial conditions. However, this 
conjecture has not yet been proven. Rather candidates of the counterexample of this conjecture have been found. 

The simplest example is the gravitational collapse of spherically symmetric dust fluid. There is an exact solution of 
Einstein equations for this case; the so-called Lemaitre-Tolman-Bondi (LTB) solution. After Eardley and Smarr have 
pointed out the occurrence of a central shell focusing naked singularity of marginally bound dust collapse , several 
theoretical efforts revealed the genericity of the shell focusing naked singularity formation in the LTB solution [|6[-^0| . 
There are several researches for more general spherically symmetric system and those revealed that non-vanishing 
pressure does not necessarily prevent the formation of the central shell focusing naked singularity pd)-p^| . 

On the other hand, there are not so much researches for non-spherical systems as the spherically symmetric cases 
since general analytic approach is impossible. Nakamura et al. and Nakamura and Sato performed numerical simu- 
lations for axisymmetric perfect fluid system pa,|6|. Their results suggest that naked singularities might be formed 



if the initial internal energy of the fluid elements is very small and initial configuration is sufficiently elongated. This 
result is consistent with the hoop conjecture proposed by Thorne, which states that black holes with horizons form 
when and only when a mass M gets compacted into a region whose circumference in every direction is C > AttM 
[fn| . Shapiro and Teukolsky also performed numerical simulations for the axisymmetric collision-less- particle system 
and showed a possibility of naked singularity formation which is also consistent with the hoop conjecture ]18|] . The 
critical behavior of the axisymmetric gravitational waves is also a candidate of the naked singularity formation [fl9| . 
An example which can be investigated analytically is the shell focusing naked singularity formation in the Szekeres 
solution |20| . Anyway, even in the non-spherical cases, there are several examples of naked singularity formation. 

There might be readers who state that numerical examples can not be regarded as candidates of naked singularity 
formation since the numerical simulation can not reveal the global structure of the space-time and further can not deal 
with infinite quantities (the numerical technique recently proposed by Hubner might be able to avoid this difficulty 
plf). However, the present authors would like to stress that if the region of sufficiently high energy, large pressure and 
large space-time curvature is visible from an observer far from the region, it should be regarded as a naked singularity. 
The region with extremely large physical quantities will be described by the quantum theory of gravity which is little 
known and hence such region is equivalent to a space-time singularity in practical sense. 

An important point of view for the naked singularity formation process has been proposed by Nakamura et al. 
[p2j| . They pointed out that the visible strong curvature region may be a candidate of strong gravitational wave 
sources. Gravitational waves generated in such a region will propagate out from there since there is no event horizon. 
From this point of view, Chiba investigated cylindrical gravitational collapse but his result is not consistent to the 
Nakamura et al's conjecture p3[ . The numerical simulations by Shapiro and Teukolsky essentially agrees with Chiba's 
result with respect to the gravitational radiation. However, we should be concerned about the numerically covered 
domain since the gravitational waves might be emitted form just the neighborhood of the naked singularity. Recently, 
present authors investigated the behavior of aspherical perturbations in the LTB space-time and showed that the 
Weyl curvature of even mode perturbations corresponding to outgoing gravitational waves diverges p3. Although 
this result was obtained by numerical integration of linearized Einstein equations, the numerical stability is guaranteed 
much better than that of the numerical simulation of full Einstein equations. Further, linearized Einstein equations 
were solved as the characteristic initial value problem and hence the numerically covered domain of the space-time 
by this analysis is much wider than the previous two numerical simulations. The results of the linear perturbation 
analysis of the LTB space-time suggest that gravitational waves will be emitted in the formation process of naked 
singularities. 

The purpose of the present paper is to re-analyze the dynamics of perturbations of the LTB space-time in the 
framework of the Newtonian approximation. In order that the singularity of the spherically symmetric space-time 
is naked, "the gravitational potential" 2M/R is smaller than unity in the neighborhood of the singularity, where 
M is the Misner-Sharp mass function and R is the areal radius. The central shell focusing naked singularity of the 
LTB space-time satisfies this condition and further the gravitational potential vanishes even at this singularity. The 
speed of the dust fluid is also much smaller than the speed of light before and at the central shell focusing naked 
singularity formation. Therefore the Newtonian approximation seems to be available, even though the space-time 
curvature is infinite at the singularity. The advantage of the Newtonian approximation scheme is that the dynamics 
of perturbations of the dust fluid and gravitational waves generated by the motion of the dust fluid are separately 
estimated; the evolution of the perturbations of the dust fluid are obtained by the Newtonian dynamics and the 
gravitational radiation is obtained by the quadrupole formula. Hence the semi-analytic estimate of the gravitational 
radiation due to the matter perturbation of the LTB space-time is possible if we adopt the Newtonian approximation. 
The results of the Newtonian perturbation analysis well agrees with the relativistic perturbation analysis by the 
present authors. This suggests that the Newtonian analysis will be a powerful tool in the analysis of some category 
of naked singularity. However, we should stress that the neighborhood of the naked singularity is not Newtonian in 
ordinary sense, because there is indefinitely strong tidal force. Further, studies of quantum effects have revealed that 
the violent quantum particle creation occurs in this space-time p5| ~ |29| . The particle creation is just a highly relativistic 
phenomenon. The Newtonian approximation scheme is available to describe the dynamics of the neighborhood of the 
naked singularity but the situation is not Newtonian. 

This paper is organized as follows. In Sec. 2, we briefly review the LTB solution. In Sec. 3, we consider the 
Newtonian approximation of a spherically symmetric dust collapse and show the relation between the Eulerian, 
Lagrangian and synchronous-comoving (SC) coordinate systems in this approximation scheme. Further in this section, 
we show the validity of the Newtonian approximation scheme even in the neighborhood of the central shell focusing 
naked singularity as long as we adopt the Eulerian coordinate system. The basic equations for the even mode 
of perturbations are presented in Sec. 4. In order to estimate the gravitational radiation generated by aspherical 
perturbations of dust fluid, we need the knowledge about mass-quadrupole moment. We show an explicit expression 
of the quadrupole moment in Sec. 5. Then we show the numerical calculation in Sec. 6 and the asymptotic behavior of 
the mass-quadrupole moment is presented in Sec. 7. Finally, Sec. 8 is devoted for summary and discussion. 
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In this article, we adopt G = c = 1 unit and basically follow the convention and notation in Ref. |3(J. The Greek 
indices mean components of tensors, while the Latin indices except for I and m represent a type of a tensor. We use 
the sub- or super-scripts I and m to denote spatial components of tensors. 

II. THE SPHERICALLY SYMMETRIC DUST COLLAPSE 

In general, the gravitational collapse of a spherical dust ball produces a shell focusing naked singularity at the 
symmetric center. This singularity can be locally or globally naked in accordance with the initial rest-mass density 
configuration. The solution of the Einstein equations describing such a situation is known as the LTB space-time. 
The line element is given by 

ds 2 = -dt 2 + {drR )} dr 2 + R 2 dQ 2 , (2.1) 
1 + f{r) 

where dfl 2 is the line element of an unit 2-sphere. The stress-energy tensor of the dust fluid is given by 

Tab = pu a u b , (2.2) 



where p is the rest- mass density and u a is the 4- velocity of the dust fluid element. In the coordinate system (2.1), the 
components of the 4-velocity is given by 

u' 1 = (1,0,0,0). (2.3) 

Einstein equations and the equation of motion for the dust fluid lead to the equations for the areal radius R and 
rest-mass density p as 

(d t R) 2 r = f(r) + ^±, (2.4) 

F(r) = 8vr ( p{d r R) t R 2 dr, (2.5) 
Jo 



where f(r) and F(r) are arbitrary functions. Eq.(2.4) might be regarded as an energy equation of the dust fluid 
element at r. From this point of view, the function f(r) corresponds to the specific energy of the dust fluid element 
and the function F(r) can be regarded as the gravitational mass function. The solution is completely fixed by choosing 
the specific energy /(r) and initial rest-mass density profile p(0,r), or equivalently the mass function F(r). The shell 
focusing naked singularity is formed only at r — 0. 

For simplicity, hereafter we will focus on the marginally bound case f(r) = 0. Since our interest is on the behavior 
of the spac e-time near the central singularity, this restriction might not lose generality of the conclusion. The solution 



of Eq.(2.4) is given by 



R=(jp) 3 [t R (r)-t]i, (2.6) 

where tn(r) is an arbitrary function which corresponds to the moment of the singularity formation. We choose the 
time of singularity formation as 

r 3/2 n[ 

tR{r) = ~\jF> (2J) 

so that R agrees with r at t = 0. 

In the spherically symmetric space-time, we can naturally introduce a "gravitational potential" <J>g defined by 

*G^|. (2.8) 

As Hayward discussed, the gravitational potential <3?g is deeply related to the formation of trapped region and hence 
gives a measure of the strength of the gravitational field |H| . To see the behavior of the gravitational potential $g 
and the velocity (dtR) r of the dust fluid element, we adopt the following initial density configuration for p 
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1 dF l | , { n \ 



where r\ and are positive constants. Since a sufficient condition for the nakedness of the centrai sheil focusing 
singularity is d 2 p(0, r) < the above initial density distribution necessarily leads to the naked singularity at the 
symmetric center r = 0. Here note that the moment of the central shell focusing singularity is t — 1 by the above 
initial density profile. The core radius of the above configuration is given by 

r G ore = n + 2 r 2- ( 2 - 10 ) 

If we set appropriate r% and T2, the space-time is globally naked singular. In Fig.l, we depict the gravitational 
potential <E>g with r\ = 2 x 1CP 2 and = 1CP 2 . We see that the gravitational potential $g is much sm aller than 



unity even at the moment of the central shell focusing naked singularity formation. Note that form Eq.( |2.4| ), \(dtR) r \ 
is also much smaller than unity. Hence in this example, the Newtonian approximation scheme seems to be available 
and in reality it is true as will be shown in the next section. 



III. NEWTONIAN APPROXIMATION 



In this section, we consider gravitational collapse of a spherically symmetric dust fluid in the framework of the 
Newtonian approximation and show that the Newtonian approximation is valid even at the moment of the central 
shell focusing naked singularity formation if the initial condition is appropriately set up as in the case of the example 
in the previous section. 



A. Eulerian Coordinate 



In the Newtonian approximation, the maximal time slicing condition and Eulerian coordinate (for example, the 
minimal distortion gauge condition) are usually adopted. The line element is expressed in the following form 

ds% = -(l + 2$ N )dT 2 + dR 2 + R 2 dtt 2 1 (3.1) 

where is Newtonian gravitational potential and we have adopted the spherical-polar coordinate system as a spatial 
coordinates. The equations for the spherically symmetric dust fluid and Newtonian gravitational potential $n are 
given by 

d T p+-^d R (R 2 pV)=0, (3.2) 
d T V + Vd R V = -d R $ N , (3.3) 
^d R {R 2 d R ^) = Airp, (3.4) 

where V is the velocity of the dust fluid element. The assumptions in the Newtonian approximation are 

\V\ < 1, and |$ N | < 1, (3.5) 

and further 

\drV\<\d R V\, |cV$n| « |9 fl $ N | and \drp\ < \d R p\. (3.6) 



B. Lagrangian Coordinate 



For the purpose to follow the motion of a dust sphere, the Lagrangian coordinate is more suitable than the Eulerian 
one. The transformation matrix between the Eulerian and Lagrangian coordinate systems is given by 

dT = dr, (3.7) 
dR = Rdr + R'dx, (3.8) 
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where regarding r and x as the independent variables, a dot means a partial derivative with respect to r while a 
prime denotes a partial derivative with respect to x. Then the line element in the Lagrangian coordinate system is 
obtained as 



dsl = - (l + 2$ N - R 2 ) dr 2 + 2RR'drdx + R dx 2 + R 2 dVt 2 . (3.9) 
Equations for the dust fluid and Newtonian gravitational potential are given by 



F(x) = 8tt / pR'R 2 dx, 
Jo 



(3.10) 



V 2 = R 2 = .f(x) + ^ : (3.11) 

*n = ^p F (*)> ( 3 - 12 ) 

where f(x) and F(x) are regarded as arbitrary functions. Since the equation for the areal radius R is the same as 
that of the LTB space-time, its solution for the marginally bound collapse f(x) = is given by the same functional 
form as Eq.(|l]), 

R ={ — ) M*)-t]3, (3.13) 

where tr(x) is an arbitrary function which determines the moment of singularity formation. 

Here we consider the Newtonian approximation of the example given in the previous section. Hence we choose the 
moment of the singularity formation as 

tr(x) = — (3.14) 

so that R agrees with x at r = 0. As for the initial density configuration, we adopt the same functional form as 
Eq.Q, 

*»•*>- £? - s { 1 + ^ (-£;) } { 1 + »» (fer) }" ■ <3 - 15) 

where xi and xi are positive constants. The above choice guarantees the regularity of all the variables before the 
singularity formation and that the central shell focusing singularity is for med a t r = 1. 

Imposing a boundary condition $n — ► for x — > oo, the solution of Eq.(3.12) is formally expressed as 

$ N = $ni(t,x) +$ N2 (t), (3.16) 

where 

$m(T,x) = J ^2 p dx, (3.17) 

^(r)^-J o ^™* = - 2 J ^ (3 - 18) 

The behavior of F/R is the same as that of 4>g as shown in Fig.l. On the other hand, the Newtonian gravitational 
pote ntial <&n of X\ = 2 x 1CP 2 and X2 = 1CP 2 is depicted in Fig. 2(a). Form these figures, we find that the condition 
( |3.5[ ) is satisfied even at the moment of the c entral shell focusing naked singularity formation. Here it is worthwhile 
to note that the right hand side of Eq.( 3.12 ) at x — diverges at the moment of the central shell focusing naked 
singularity formation, 

*n — > ~m %~ 1/3 for x — ► at t = 1, (3.19) 
27 V ' 

where 
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1 d 2 T R (x) 

2 dx 2 



(3.20) 



However since the power index of x is larger than —1, <E>n itself is finite at x = even at the moment of the central 
shell focusing singularity formation r = 1. 

In order that the Newtonian approximation is successful, temporal derivatives of all the quantities should always 
smaller than the radial derivatives of those. Here we shall focus on the neighborhood of the central shell focusing 
naked singularity only. For this purpose, we introduce a new variable w defined by 



w = St l ' 2 i 



(3.21) 



where St = (1 — t)/^^). Then we consider a limit r — > 1 with fixed w. It should be noted that x also goes to zero 
by this limiting procedure. The mass function F, rest-mass density p and areal radius R behave as 



F 

P 
R 



-w 3 St 3 ' 2 , 



2tt 



(3.22) 

(3.23) 
(3.24) 



All these variables are proportional to the power of St and the coefficients of those are functions of w. It is easy to 
see that their derivatives with respect to r or x also take the same functional structure with respect to St and w. 
Thus the ir-dependent part <&ni of the Newtonian gravitational potential will also behave in the manner 



$ni — ► n1 (w)St\ 



(3.25) 



where i is a constant and </>n i s a f unction of w. Substituting the above equation into Eq.(3.12) and using the 
asymptotic behavior (3.22) and (3.24), we obtain 



dcj) m (w) 
dw 



2w(3 + 7w 2 



27^(1 + ^)5/3 



St- 1 '*. 



(3.26) 



In order that the depen denc e of both sides in the above equation on St agrees with each other, i should be equal to 
1/3. Integration of Eq.( 3.26 ) leads to 



(f>m(w) 



1 



9^(1 + ^)2/; 



{7(l + w 2 )-9(l + w 2 ) 2/3 + 2} 



(3.27) 



In order to see the asymptotic dependence of 3>n2{t) on r, we differentiate Eq.(3.18) to obtain 



X2 



F' F , 
R~^R dX - 



(3.28) 



The integrand in the right hand side of the above equation behaves near the origin at the moment of the central shell 
focusing singularity formation as 



F' F 
R?\R 



-7/3 



9t 



5/3 
R(l) 



for 



at r = 1. 



(3.29) 



Therefore, the integral in Eq.(3.28) does not have finite value at r = 1. Since as shown in the above, this divergence 
comes from the irregularity of the integrand at the origin x — 0, we shall estimate the contribution near the origin to 
the integral in Eq.(3.28). We again consider the limit of St — > with fixed w and obtain 



F' IF 
IPVR 



dx 



9r 



5/3 
R(l) 



■St- 2 ' 3 



'jdw 



(l + W 2 ) 5 /3 



3t 



5/3 
Ml) 



■St- 2 ' 3 . 



(3.30) 



Substituting the above equation into Eq.(3.28) and integrating it with respect to r, we obtain 



G 



St 1 / 3 

$N2 — » ^ 7r + $N 2 (l)- (3.31) 
T fl(l) 

As a result, in the limit of <5t — > with fixed to, $n is expressed in the form 

$ n — > 2/ f + 7w2 ; St 1 / 3 + $ N2 (1). (3.32) 

9^(1 + ^)2/3 

We depict numerically obtained $n together with the above asymptotic form in Fig. 2(b). The asymptotic estimate 
agrees with the numerical result quite well. 

Now we have known that in the limit of St — > with fixed to, all the variables behave as 

Z(t,x) — > z(w)St 3 + constant, (3.33) 

where z(w) is some function of to and j is a constant. The derivatives of Z with respect to T and R are expressed by 
using its derivatives with respect to r and x as 

(d T Z) R = (d T T) R Z + {d T x) R Z', (3.34) 

(d R Z) T = {8 r t) t Z + (d R x) T Z'. (3.35) 



From Eqs.(^7J) and fl3.8|) , we find 



R 1 

(9 T T) fl = 1, (d R T) T = 0, (d T x) R = -— and (9ija;)T = -~ r (3.36) 



Then the following relation is derived, 



(d R z) T z 

Inserting Eq.( 3.33| ) into the above equation, we obtain 



{d T z) R _ Rl z__ k (33?) 



-l 



{9tZ)r • ' ll w{1 + 7w 2)- jz{ 3 + 7w*)(¥-\ } StK (3.38) 



(d R Z) T 3r i/3 )(1 + u; 2)i/3 |2 J \dwJ 

The above equation means that in the limit of St — > with fixed to, the following inequality holds 

\{d T Z) R \ « \(d R Z) T \. (3.39) 

Therefore, the order counting of the Newtonian approximation is guaranteed even in the neighborhood of the central 
naked singularity. 

Here it is worthy to notice that in the limit of St — > with fixed to, \Z\ is much larger than \Z'\ because of 



Z 1/1 dw 

-to — jz- 



Z' t R (d V 2 dz 



8t- 1/2 — ► oo. (3.40) 



Hence the vicinity of the central shell focusing naked singularity is not Newtonian situation in ordinary sense. In 
reality, the present authors have shown that the violent quantum particle creation occurs in such a situation p8| , p9| . 
This means that the formation of the central shell focusing singularity is a highly relativistic phenomenon even though 
it is well described by the Newtonian approximation in the Eulerian coordinate system. 



C. Synchronous-Comoving Coordinate 

In the LTB solution, the SC coordinate system is adopted. We should note that the SC coordinate system is 
different from the Lagrangian one. We consider the following coordinate transformation 
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dt 



1 + 



1 



{d T R)%\ dr - (d T R) x (d x R) T dx, 



dr = dx. 



(3.41) 
(3.42) 



Assuming Eq.( |3.5| ), derivatives of the areal radius R with respect to the Lagrangian time coordinate r and the radial 
coordinate x are written as 



R = {d r t) x (d t R) r + (d T r) x (d r R)t = {1 + $ N - {d T R)l} (d t R) x ~ (d x R) T , 
R' = (d x t) T {d t R) r + {d x r) T (d r R) t = {d r R) t - {d T R) x {d x R) T {d t R) r ~ {d r R) t . 



(3.43) 
(3.44) 



Hence the line element in the new coordinate system t, r up to the lowest order in the sense of Eq^SJj) is written as 

-dt 2 + (d r R) 2 t dr 2 + R 2 dQ 2 . (3.45) 



dSgQ 



From the above equation, we find that the coordinate transformation ( 3.41 ) and ( 3.42| ) leads to the SC coordinate 
system. Especially, the above line element is completely the same as the relativistic one in the marginally bound case. 
Further we obtain the equations of the lowest order in the SC coordinate system are given in completely the same 
form as that in the relativistic one: 



{d t Rf r = f(r) + 



F(r) 
R ' 



(3.46) 
(3.47) 



F{r) = 8tt f p(d r R) t R 2 dr. 
Jo 

It is worthwhile to note that in the SC coordinate system, the Newtonian gravitational potential $n does not appear 



IV. BASIC EQUATIONS OF EVEN MODE PERTURBATIONS 

We consider non-spherical linear perturbations in the system of the spherically symmetric dust ball described in 
the previous section. First, we consider perturbations in the Eulerian coordinate system. The line element is written 
as 

ds% = - (1 + 2$ N + 2£$ N ) dT 2 + dR 2 + R 2 dn 2 , (4.1) 



wher e 5$n is a perturbation of the Newtonian gravitational potential. Using the transformation matrix (3.7) and 
(3.8), we obtain the perturbed line element in the background Lagrangian coordinate system as 



dsl = -(l + 2$ N + 25$ N - R 2 ) dr 2 + 2RR'drdx + R' 2 dx 2 + R 2 dtf. (4.2) 



Hereafter we discuss the behavior of perturbations in this coordinate system. The density p and 4-velocity u M are 
written in the form 

p = p(l + 6 p ), (4.3) 

= u" + Su". (4.4) 

By definition of the Lagrangian coordinate system, the components of the background 4-velocity is given by 

(«") = (tZ , 0,0,0). (4.5) 

From the normalization of the 4-velocity, we find 

5u° = -<5$ N + RR'Su 1 . (4.6) 

The order-counting with respect to the expansion parameter e of the Newtonian approximation is given by 

Su" = 0(e 2 ), 6u e = 0{e), S p = 0(e°) and 8<S> N = 0(e 2 ). (4.7) 

Then the equations for the perturbations are given by 
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1 



d T 5 p + -=-^di (p^/j5u e ) = 0, 

0. 



d T Sue 



N = 



de {V^l lm d m 5<S> N ) - inpSp = 0, 



where 



^=R'R 2 sm9, 



(4.8) 
(4.9) 
(4.10) 

(4.11) 



and r ) lm is a contravariant component of the background 3-metric. 

Here we focus on the axisymmetric even mode of perturbations. Hence the perturbations are expressed in the form 



S p =^A p (i)(r, x)P t (cos8), 
l 

5<f> N = ^ A m (T,x)Pi(cos6), 
i 

Sui = ^2u x{1) (t,x)P 1 (cos8), 
i 

6u 2 = Y u e(i)( T ,x)^Pi(cos6) 1 
i 

5u 3 = 0. 



From Eqs.(f4.8|), (|4.9|) and (|4.10| ), we obtain 

A p(0 + 



F' \R' 2 



R 2 

U x{ i) + A' $(z) 



U$(i) + A$(j) = 0, 



R'R 2 



= 0, 



(4.12) 
(4.13) 
(4.14) 

(4.15) 
(4.16) 

(4.17) 

(4.18) 
(4.19) 

(4.20) 



Comparing the basic equations for the relativistic perturbations of I — 2 in Ref. J24| to the above equations, we find 
the correspondence between the Newtonian and relativistic variables; A$( 2 ) = — 2K, U x n) = — V\ and f/g(2) = — Vi. 



V. MASS-QUADRUPOLE FORMULA 

Hereafter we focus on the quadrupole mode I = 2 and therefore omit the subscript (I) to specify the multi-pole 
component of the perturbation variable. The mass-quadrupole moment Qi m is given by 



Q em = J p \XtX m - ^R 2 S em ^j d 3 X = ^Q(T)diag[-l, -1, 2] 



where 



Q(T) 



pApR dR. 



For a function g(r, x) with sufficiently rapid fall off for x — * oo, we find that 

/>oo />OC /• 

■ / gdR = / (d T g) R dR = / 
Jo Jo Jo 



d_ 

dT 



R' 



r°° / . \ r . -i oo r°° 

J [R'g + R'gj dx - \Rg\ ^ = J d T (R'g)dx 



(5.1) 
(5.2) 

(5.3) 
(5.4) 
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Using the above formula, we obtain 

jmn /'CO om i poo om 

I -k— (p& P R'R 4 ) dx=— F'-^- — (A p i? 2 ) dx. (5.5) 
dT m J dr m VP 9 ' 8vr J 9r m v p ; v ; 

The power Lqw carried by the gravitational radiation at the future null infinity T + R — ► oo is given by 

32tt 2 f rf 3 Q( M ) V 

iGW = l75-(^^J ' (5 - 6) 

where u = T — R is the retarded time. The Weyl scalar $4 carried by outgoing gravitational waves at the future null 
infinity is estimated as 

*4 = -C abcd n a m b n c fh d = sin 2 9, (5.7) 

5 our 

where C a bcd is the Weyl tensor, and n a and m a are two of the null tetrad basis whose components in the spherical 
polar coordinate are given by 



1 
1 

7t v 



K) = — (1,-1,0,0), (5.8) 
(m M ) = -^=(0, 0, R, -iRsmd). (5.9) 



It should be noted that the power I/gw is proportional to the square of the 3rd-order derivative of Q(u) while the 
Weyl scalar ^4 is proportional to the 4th-order derivatives of Q(u). 

VI. NUMERICAL SIMULATION 

In this section, we explain the procedure to numerically integrate the system of partial differential equations for 
the quadrupole perturbation variables and show the time evolution of the mass-quadrupole moment. 

A. Basic Equation 

We assume that all the perturbation variables are regular before the central shell focusing singularity formation 
and hence are written in the form of the Taylor series as 



A p = x 2 A; = x 2 Y,^ P (n)X 2n , (6-1) 

00 

U x = xU:=xJ2u* x (n)X 2n , (6-2) 

CO 

Ue^x'U^x^U^x 2 -, (6.3) 

00 

A* = ,t 2 A^ = x 2 A* (n) x 2 ". (6.4) 

71=0 

In the numerical calculation, we focu s on the above asterisked variables rather than the original variables. 



To obtain the solution of Eq. ( 1.2C ) , we introduce a variable Q defined by 



A* = ^ 3 . (6.5) 



Substituting the above form of A$ into Eq.(4.20), we obtain the equation for Q as 
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Integrating the above equation, we obtain 



■A p (r, xx), 



# 3 (t,zi) pv 

where we have chosen the integration constant so that Q' is finite for x — ► cxd. Further integration leads to 

- / dx l R'{T,x 1 )R 4 {T,x 1 ) / dx 2 v 2 \ Ap(r ; : 

2 JO Jx R\T,X 2 ) 



(6.7) 



1 r x r°° f'(t ) 

Q{r,x) = -- dxx&faxj&frx!) / dx 2 p3 / 2 \ A p (t,x 2 ) 

2 Jo i? d (r,x 2 ) 



,£2) 



1 r dx n**) 



R 3 {t,x 2 ) 



A p (t, cc 2 ) / dxii?'(r,a;i)i? 4 (T,a;i) 



where 



= - I Q{i? 5 -4(r,x)+ a; 7 S(T,x)} 



A(r,a;)= / ttoi J^^ 1 ] . A p (r, a?i ) , 



1 f x 
x Jo 

and we have chosen the integration constant so that Q vanishes at the origin, x = 0. Hence A J is written as 



1 (R 



x \ ' 



(6.8) 

(6.9) 
(6.10) 

(6.11) 



For the numerical calculation, we rewrite the basic equations (4.17) — (4.19) for the p ert urba tion variables. Differ- 
entiating the conservation law ( 4.17 ) with respect to time r and then using Eqs.( 4.17 )— ( 4.19 ), we obtain a second 
order differential equation with respect to r for A* as 



dv; 

dr 

sa; 

dr 



A, 



R' p 2R'R 2 



F' 2 
A* 



9 xR' 2 V R 



R \ TT , 12 R R'\ 

-1 u *-rAr-r>) u ° 



R? 



xF'R V R'R 2 



F> V 



= v* 



where 



D| = Va£)' = — r//>> 

a; 10 



The equations for U* and t/g are written as 



dr 

mi 

dr 



- r>* 



-A* 



(6.12) 

(6.13) 
(6.14) 

(6.15) 

(6.16) 
(6.17) 



Eqs.(3.11) — (6.17) constitute a closed system of equations, which is solved numerically. 
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B. Background and Initial Condition 



The background solution is given by (3.13), (3.14) and (3.15). The parameters x\ and xi in p(0,x) are chosen 
to be x\ = 2 x 10~ 2 and X2 = 1CP 2 , respectively. Then the core radius x COYC of the background rest-mass density 
distribution becomes 



x x + -x 2 = 1.25 x 10~ 2 . 



(6.18) 



We start the numerical integration at r = 0. The initial data of the density perturbation A* is set up in the form 

2 ~l 



A* p (0,x) =exp- 



1 / x 



2 V (TpX c 



where a p is a positive constant smaller than unity. We set a p = 0.5. 
From Eqs.fl4.18|) and ( [4.19[ ), we find 

d T (U x - U' e ) = 0. 

Hence U x is written in the form 

U x = U' e + C x {x), 

where C x is an arbitrary function. As for the initial data, we set Ug = and 



C x = c x x exp < -- 



1 / x 



2 V *^ C 



(6.19) 

(6.20) 
(6.21) 

(6.22) 



where a x is a positive constant smaller than unity and c x is also a constant. From Eq.(4.17), the initial data for V* 
should be 



1 



v p *(o,x) = --^(F'c x y 



(6.23) 



C. Numerical Scheme 



The temporal integration of Eqs.(6.12)— (6.17) i s pe rform ed by the second order Runge-Kutta method. In order to 
obtain A|, and DJ, we numerically i ntegr ate E qs.(|6.9| ) and (3.10) also by the second order Runge-Kutta method and 
then substitute the results into Eqs.( 6.11 ) and (j6.15| ). 

The number of spatial grid points is 5 x 10 4 and the outermost grid point corresponds to x — 10x c . Since we 
numerically integ rate the second order differential equation with respect to r for the d ensity perturbation A*, the 
conservation law (4.17) is not trivially satisfied by the numerical error. Hence we use Eq.(4.17) as a check of accuracy 
of the numerical integration. We define the error function £ as 



£ 



V* 



zE 12 



U* 



( xF 1 



6 

i? 2 



\v* 



u* 



cR' 2 



v: 



x 2 F' 



xF' 

Iv 2 



6 

R 2 



(6.24) 



We perform the numerical integration until 1 — t = 2 x 10~ 4 . The function £ in x < x c is less than 10 3 . Hence 
the numerical integration is almost consistently performed near the origin and hence we can observe the asymptotic 
behavior of the contribution near the central shell focusing singularity to the mass-quadrupole moment in the limit 
ofr->l. 



D. Behavior of the Mass-Quadrupole Moment 

In Fig3, we show d 2 Q/dT 2 , d 3 Q/dT 3 and d i Q/dT i . The numerical integration was performed from 1 — r = 1 to 
1 - r = 2 x 10~ 4 . From this fi gure, we find that those variables do not show power-law behavior with respect to 
1 — r but rather their behavior is oscillatory and growing. At this stage, the mass-quadrupole moment is determined 
mainly by the motion of dust fluid in the outside region and the contribution in the neighborhood of the central shell 
focusing singularity is still negligible. However, the contribution near the origin will be dominant in the limit r — > 1 
as will be shown in the next section. 
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VII. ASYMPTOTIC ANALYSIS OF THE PERTURBATIONS 



In order to get information of the asymptotic behavior of the mass-quadrupole moment, we should carefully examine 
the asym ptotic behavior of the perturbation variables near the origin. For this purpose, we introduce w defined in 
Eq.( |3.2l' ) and then consider the limit r — » 1 with fixed w. Since all the background variables appearing in the 
equations of th e pertu rbatio ns ar e proportional to the power of St and the coefficients of those are the functions of 
w as Eqs. ( |3.22] ) — ( 3.24 ) and ( 3.32 ), we expect that the perturbation variables also behave in the same manner as the 
background variables and hence we assume 



v; = s;(w)6t->, 
i 



U* x = -d x {x 2 U* e ) = 7$) + 2u* g ) Sr ". 

u e = t r{i) u K w ) St ~ 9 ') 



du* D 



^ - t r{() 5 1{ w ) 5t r > 



A* — 



(7.1) 
(7.2) 

(7.3) 
(7.4) 



where we have assumed that the contribution of C x to U* is negligible in the limit of r — > 1. In Fig. 4, we depict 
the numerical results for the asterisked variables at the origin as functions of 1 — r in the case of C x = 0. From this 
figure, we can easily find that the hatted variables show the power-low behavior and grow monotonically. The power 
indices are estimated as 



p = 1.70, q = 0.38 and 



1.37. 



(7.5) 



Due to this growing behavior, in reality, the contribution of C x to U* becomes negligible in the limit r — ► 1. Although 
we do not show the case of C x 0, we have confirmed that the behavior of the asterisked variables in the case of 
C x is almost identical to the case of C x — in the limit r — > 1 . 

Next we depict the following normalized perturbation variables with respect to w for x < 0.05a; core at various time 
steps for the time interval 2.0 x 10 -4 < 1 — r < 6.6 x 10~ 3 in Figs. 5, 6 and 7; 



At=A;(r,i)/A;(r,0), 
U$ = U2(t,x)/US(t,0), 
Ai = A%(r,x)/A%(r,0). 



(7.6) 
(7.7) 
(7.8) 



Note that A*(t, 0) at 1 — r = 2 x 10~ 4 becomes 381 times larger than that at 1 — r = 6.6 x 10~ 3 . As expected, 
t he sp atial configuration of these variables with respect to w are almost time independent. Hence the assumptions 
([7.l|)— ( |7^) are justified by the numerical calculation. 

Now, by virtue of our knowledge about the asymptotic forms (7.1 ) — (|7.4|), more rigoro us di scu ssion about the 



evolution of the mass-quadrupole moment is possible. Substituting Eqs. Q7.ip — (7.4) into Eqs. (4. 17) — (4.20), and using 
the asymptotic behavior of the background variables ( [3. 22 ) — (3.24), we obtain 



' d6 *o 

w—^ + 2p5) 
. aw h 



--p-i 



18 d 

w A dw 



w 2 (l + w 2 ) 2/3 d 
(3 + 7w 2 ) 2 ~d*w 



(w 2 u* 9 ) 



Ylul 



w 2 (l + w 2 ) A / 3 



Sr-1- 7 / 3 = 0, 



(7.9) 



w _^£ + 2qu* g ) St-"- 1 + 28%8T- r = 0, 
dw ' ' 



and 



_d_ 

dw 



w 2 (l+w 2 ) 5 / 3 d 
3 + 7w 2 dw 



(w 2 8l) 



2w 2 (3 + 7w 2 ) 

3 (l + w 2)l/3 ' 



5r -r+5/3 _ _ w 4 S * St -P+2 = Q 



Since the power of St should be balanced in each equation, we obtain 

4 . 1 



P 



and 



r = p - 



(7.10) 



(7.11) 



(7-12) 



Eqs. (7.9)— (7.11) constitute a closed system of ordinary differential equations. In order to obtain a solution of the 
above equations, we also need numerical integration. However, here we are going to search for regular and gentle 
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solutions of the ordinary differential equations. This is much easier than the previous numerical integration of the 
partial differential equations of which solutions show singular behavior. Further high numerical accuracy is guaranteed 
and therefore the following analysis might be called "semi-analytic" . 

By an appropriate manipulation |32| , we obtain a single decoupled equation for u* g as 



where y = w 2 and 



dy 4 



c 3 - 



dy 3 



d 2 ut dut 



Co 



= -2{-24(9 + 26y + 21y 2 ) 

+ 3q 2 (63 + 414y + 1016y 2 + 1106y 3 + 441y 4 ) 

+ g(252 + 1323y + 3149y 2 + 3465y 3 + 1323y 4 )} 

x {92/ 3 (l + y) 3 (3 + 7y) 2 }~\ 



(7.13) 



(7.14) 



ci = {378 + 2196y + 4758y 2 + 2006y 3 - 4424y 4 - 3234y 5 
+ 3q 2 (l + y) 2 (189 + 1035y + 1911y 2 + 1225y 3 ) 
+ g(1323 + 7893y + 18966y 2 + 21202y 3 + 9247y 4 + 441y 5 )} 
x {18y 3 (l + y) 3 (3 + 7y) 2 }-\ 



(7.15) 



c 2 = {1269 + 7731y + 18453y 2 + 19565y 3 + 7350y 4 + 9q 2 (3 + Wy + 7y 2 ) 2 
+ 6q(153 + 993y + 2387y 2 + 2527y 3 + 980y 4 )} 
x {9y 2 (l + y) 2 (3 + 7y) 2 }- 1 , 



(7.16) 



c 3 = {159 + 506y + 427y 2 + 12g(3 + lOy + 7y 2 )} 
x {6y(l + y)(3 + 7y)}- 1 . 



(7.17) 



Giving an appropriate boundary condition, we can numerically solve Eqs. (7.9) — (7.11) as a kind of the eigen value 
problem to obtain q and the solution for u* g . 

The boundary condition at y = is given by the following procedure. First, we assume that 5*, u g and <5J are C°° 
functions. By this assumption, the perturbation variables are written in the form of the Taylor series as, 



(m)D 



* = E 



u( m )y 1 



m=0 



51= E A 



{m)V 



m— 



Then substituting the above expressions into Eqs. ( 4. 17 )— ( 4.20D , we obtain 



<5(o) ( 3 + 9 ) - 32u (0 ) + 14u (1) = 0, 



S W [ 3 + '-<' 



''(2) 



10 



18848 



1568 

7 

9 

1388 



*(0) 



Ho) 



104u (1) + 36u (2) = 0, 

(3) = 0, 

(q + m)u (m) + A (m) = 0, 



um — 208u( 2 ) + 66u 



(7.18) 
(7.19) 
(7.20) 

(7.21) 
(7.22) 
(7.23) 
(7.24) 



-28A 



(o) 



21A ( i) - 5(o) = 0, 



(7.25) 



14 



148A (0 ) - 3(46A (1) - 54A (2) + 5 (1) 



620 

-^T A (o) 



119 
~L8~ 



A (1) - 4A 



11 



(2) 



H3) 



18 



'(2) 



0, 
0, 



where m is non-negative integer. From the above equations, we obtain 
4(-24 + 28q + 21<7 2 ^ 



«(2) 



it 



(3) 



21(2 + 7 9 + 3g 2 ) " M(0) ' 
2(5248 + 4822g + 3765g 2 + 2646g 3 -I- 567q 4 ) 

567(2 + 7q + 3g 2 )(12 + 13q + 3q 2 ) M(0) ' 
4(2174336 + 3606320g + 2343504g 2 + 889242g 3 + 272511<7 4 + 62370g 5 



6237g 6 ) 



18711(2 + 7q + 3 9 2 )(12 + 13q + 3q 2 )(28 + 19q + 3q 2 ) 



Ho)- 



(7.26) 
(7.27) 

(7.28) 
(7.29) 
(7.30) 



Since Eq.( |7.13 ) is linear, the value of U(q) has no special meaning and hence we set U( ) = 1, for simplicity. Then the 
boundary condition at y = for Eq.(7.12) is uniquely determined as 



«fllw=o = «(o) = l, 



y=0 



dy 2 



2u 



(2)- 



and 



d 3 i 



dy 3 



6m 



(3)- 



(7.31) 



We numerically integrate Eq.( 7.13 ) outward from y = by 4th-order Runge-Kutta method. The behavior of Uq 
depends on the value of q. 

In order to set up the boundary condition at the outer numerical boundary, we consider the behavior of Eq. ( 7. IS ) 
in the limit of y — > oo, 



1 



dh 



1 



^ + % (129 + 61) ^# + V (3 ^ 

Therefore, u% behaves as 



, d 2 u* B 1 , . dut 6 



0. 



IT 



Const, x y 



for 



(7.32) 



(7.33) 



Inserting the above equation into Eq.(7.32), we obtain a 4th-order algebraic equation for k. The solutions of this 
equation are given by 



9 4 

k = -q, -(q+ 1), - - and - 



(7.34) 



In order to pick up a solution from the above four independent asymptotic solutions, we consider Ug in the limit of 

V -> oo, 



Un 



Const, x y k ST- q = Const, x x 2k 5T~( k+q) . 



(7.35) 



Here we impose a condition that Ug is non-zero and finite for < x < e at St = 0, where e is a positive infinitesimal 
number. This condition leads to 



k = —q. 

Hence by varying q, we search for the solution which behaves 

Ug ► Const, x y~ q for y > oo. 

The numerical calculation reveals that the above behavior is realized when 

q = 0.3672 . 



From Eq.( 7.12 ), we obtain 



p = 1.701 and 



1.367 



(7.36) 
(7.37) 
(7.38) 
(7.39) 



The above values agree with the numerical results (|7.5| ) quite well. We depict the data for 8*/8{q), u* g and <$J/A( ) in 

Figs. 5, 6 and 7 together with the variabl es At , Ug a nd A ^ obtained by the numerical calculation of Eqs.( 3T^ )— ( |6~T7 ), 
where 8t Q \ and A( ) are defined by Eqs.( 7.1q ) and ( 7.20| ), respectively. We see the quite excellent agreement. 
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Now we examine the mass-quadrupole moment Q(T) and its time-derivatives d m Q/dT m . In order to se e th e 
contributio n of th e cen tral si ngu larity to d m Q/dT m , we consider the integrand in the right hand side of Eq. (|5.5|) . 
Using Eqs.(|3.22D, (|3.24Q and (|7\l]), we obtain 



F'ApR 2 



From the above equation, we obtain 



(7.40) 



r(m) 



(r, x) 



fxm q2 — m 

U „2~\ , f 4/ 3 - m A^ 13 / 3 -P- m „, 26 /3- 2 P- 2m 

3~ T R(1) 0T " 



dr 



(F'A P R 2 ) 



x w 



dw 



)m 
| w 2(p-4/3) (1+u; 2 ) 4/3 (5 * H | 



(7.41) 



We consider the integral of I*" 1 ) from x — to x — x Q to see the contribution of the central shell focusing 
naked singularity to the time derivatives of the mass-quadrupole moment. Here we take a limit St — > with fixed 
w = XqSt^ 1 / 2 and then consider the limit w a — > oo. As a result, we obtain 



= St 1 / 2 


rw 

/ j(™ 




Jo 


2 2-m 


4/3-m, 


3 


T fl(l) 



/((' 



26/3-2p-2m 



{«,2(P-*/3)( 1 + wY^MW 



(7.42) 



The above equation and Eq.(7.39) show that the contribution of the central singularity to d m Q/dT m diverges for 
r — > 1 if and only if m is larger than or equal to 4. This result and the quadrupole formula imply that the metric 
perturbation corresponding to the gravitational radiation and its first order temporal derivative is finite but the second 
order temporal derivative diverges. Hence the power Lgw of the gravitational radiation is finite but the curvature ^4 
carried by the gravitational waves from the central naked singularity diverges; see Eqs.(|J) and Q. This agrees 
with the relativistic perturbation analysis. Further, we find that in the limit of r — ► 1, 



* 4 = - 



3ird 4 Q 



S T 5 / e - p oc (1 - r)-° 



867 



5 dT 4 

This result is also consistent with the relativistic perturbation analysis |24 



(7.43) 



VIII. SUMMARY AND CONCLUDING REMARKS 



We analyzed the even mode perturbations of I — 2 in the spherically symmetric dust collapse in the framework 
of the Newtonian approximation and estimated the gravitational radiation generated by these perturbations by the 
quadrupole formula. Since we treat separately the dynamics of the matter perturbations and the gravitational waves 
in the wave zone, we can estimate the asymptotic behavior semi-analytically, where "semi-analytically" means that 
we know it by solving the gentle ordinary differential equations. This is the great advantage of the Newtonian 
approximation. 

As a result, we found that the power carried by the gravitational waves from the neighborhood of the naked 
singularity at the symmetric center is finite. However, the space-time curvature associated to the gravitational waves 
becomes infinite in accordance with the power law. This result is consistent with recently performed relativistic 
perturbation analysis by Iguchi et al. ]24f| . Furthermore, the power index obtained by the Newtonian analysis also 
agrees with the relativistic perturbation analysis quite well. 

The agreement between the results of the Newtonian and relativistic analyses implies that the perturbations them- 
selves are always confined within the range to which the Newtonian approximation is applica ble. H ere we will focus 
on the metric perturbation, AJ. Since the asymptotic solution of A J has the same form as Eq.( 3.33] ), we immediately 
find that in the limit of St — > with fixed w, 



1G 



St 1/6 , (8.1) 



and hence the assumption |9rAJ| <C |c?i{AJ| of the Newtonian approximation is valid in the Eulerian coordinate 
system. We can also see 2nd-order derivatives. In the same limit, we find 

"W afr and wj Kfr ■ (8 ' 2) 

From the above equations, we obtain 



g|A| 



oc 5r 1/3 . (8.3) 



The above equation means that in the limit of St — > with fixed ty, the following inequality is also satisfied, 

\d 2 T A* | « |c^A* |. (8.4) 

The above inequality implies that the wave equation for the metric perturbation A$ is well approximated by the 
Poisson-type equation if we adopt the Eulerian coordinate system. 

However as mentioned in Sec. 2, the gravitational collapse producing the shell focusing globally naked singularity is 
not Newtonian in ordinary sense. The same is true for the perturbation variables because |A$/A^| 3> 1 in the limit 
of St — > with fixed w. Even though the Newtonian approximation is valid for the Eulerian coordinate system, the 
Newtonian order counting breaks down if we adopt the Lagrangian coordinate system as the spatial coordinates. 
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FIG. 1. The spatial configurations of the gravitational potential <J>g at t = 0, 0.5 and at the moment of the central shell 
focusing singularity formation t = 1 are plotted. In all the cases, the gravitational potential <£>g is everywhere much smaller 
than unity. 
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FIG. 2. (a) The spatial configurations of the numerically obtained Newtonian gravitational potential $n are plotted as a 
function of x. The solid line represents $n at r = 0, while the dashed line corresponds to the case at r = 0.95. The dotted 
fine represents the case at the central shell focusing singularity formation t = 1. Although the regularity of $n at the origin 
breaks down at t = 1, |$n| is everywhere much smaller than unity, (b) The asymptotic behavior of the Newtonian gravitational 
potential <!>n is plotted as a function of w for r = 0.999 (the t op), r = 0.99945 (the middle) and r = 0.9999 (the lowest). The 
solid lines represent <&n obtained by the asymptotic analysis (3.32) while dashed lines is the numerically obtained ones. In all 



the cases, numerically obtained values well agree with those obtained by the asymptotic analysis. 
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FIG. 3. d 2 Q/dT 2 , d 3 Q/dT 3 and d 4 Q/dT 4 are plotted as a function of 1 — r in the logarithmic scales. The absolute values 
of these quantities show oscillatory behavior and grow. At this stage, the contribution near the origin is not dominant. 
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FIG. 4. The time evolutions of the asterisked perturbation variables, A*, J7*, Ug and AJ, at the origin are plotted in the 
logarithmic scales. All these quantities monotonically grow asymptotically in accordance with the power low. 
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FIG. 5. Spatial configurations of the normalized density perturbation Aj, are plotted as a function of w for x < 0.05a; cor c at 
various time step form 1 — r = 6.6 x 1CP 3 to 1 — r = 2.0 x 1CP 4 by solid lines. The density perturbation becomes about 381 
times during this time interval. On the other hand, 8* (w) /<5( ) obtained by the asymptotic estimate is represented by a dashed 
line. It is very difficult to distinguish these lines from each other. 
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FIG. 6. Spatial configurations of the normalized tangential velocity perturbation \j\ are plotted as a function of w for 
x < 0.05a; corc at various time step form 1 — r = 6.6 x 10~ 3 to 1 — r = 2.0 x 10~ 4 by solid lines. On the other hand, Ug(w) 
obtained by the asymptotic estimate is represented by a dashed line. We can see that Ul asymptotically approaches to u* e . 
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FIG. 7. Spatial configurations of the normalized potential perturbation AJ are plotted as a function of w for x < 0.05x CO re 
at various time step form 1 — r = 6.6 x 1CP 3 to 1 — r = 2.0 x 10~ 4 by solid lines. On the other hand, 5|,(w)/A(o) obtained by 
the asymptotic estimate is represented by a dashed line. The normalized potential perturbation A^ asymptotically approaches 
to 5p/5(o). 
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